# -*- coding: utf-8 -*-
"""
文件创建时间: 2025-11-05
文件功能描述: 稳态求解日志分析工具
- 解析求解器日志文件
- 提取残差和关键数据
- 生成收敛趋势可视化图表
- 自动查找最新日志文件进行分析
"""

import re
import numpy as np
import matplotlib.pyplot as plt
from pathlib import Path
from typing import Dict, List, Tuple
import pandas as pd
import os
import sys

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'Microsoft YaHei']  
plt.rcParams['axes.unicode_minus'] = False

# 添加自动查找最新日志文件的功能
def find_latest_log():
    """查找最新的日志文件"""
    log_dirs = [
        Path(r"e:\mycode\subchannel\TwoFluidSCA\results"),
        Path(r"e:\mycode\subchannel\TwoFluidSCA\src\results"),
    ]
    
    latest_file = None
    latest_time = 0
    
    for log_dir in log_dirs:
        if log_dir.exists():
            for log_file in log_dir.glob("*.log"):
                mtime = os.path.getmtime(log_file)
                if mtime > latest_time:
                    latest_time = mtime
                    latest_file = log_file
    
    return latest_file

class ConvergenceAnalyzer:
    """收敛分析器：解析日志并可视化收敛过程"""
    
    def __init__(self, log_file_path: str):
        """
        初始化分析器
        
        Args:
            log_file_path: 日志文件路径
        """
        self.log_file = Path(log_file_path)
        self.data = {
            'iterations': [],
            'res_vapor_axial': [],
            'res_liquid_axial': [],
            'res_vapor_lateral': [],
            'res_liquid_lateral': [],
            'res_mass_vapor': [],
            'res_mass_liquid': [],
            'res_pressure': [],
            'res_void_fraction': [],
            'diag_dominance_vapor': [],
            'diag_dominance_liquid': [],
            'relaxation_pressure': [],
            'relaxation_velocity_g': [],
            'relaxation_velocity_l': [],
        }
        
    def parse_log(self):
        """解析日志文件，提取关键数据"""
        print(f"开始解析日志文件: {self.log_file}")
        
        if not self.log_file.exists():
            raise FileNotFoundError(f"日志文件不存在: {self.log_file}")
        
        with open(self.log_file, 'r', encoding='utf-8') as f:
            content = f.read()
        
        # 正则表达式模式
        patterns = {
            'iteration': r'第\s+(\d+)\s+次迭代:',
            'res_vapor_axial': r'气相轴向动量残差:\s+([\d.e+-]+)',
            'res_liquid_axial': r'液相轴向动量残差:\s+([\d.e+-]+)',
            'res_vapor_lateral': r'气相横向动量残差:\s+([\d.e+-]+)',
            'res_liquid_lateral': r'液相横向动量残差:\s+([\d.e+-]+)',
            'res_mass_vapor': r'气相质量守恒残差:\s+([\d.e+-]+)',
            'res_mass_liquid': r'液相质量守恒残差:\s+([\d.e+-]+)',
            'res_pressure': r'压力残差:\s+([\d.e+-]+)',
            'res_void_fraction': r'空泡份额残差:\s+([\d.e+-]+)',
            'diag_dominance': r'对角占优度:\s+([\d.e+-]+)',
            'relaxation_pressure': r'压力松弛因子至([\d.]+)',
            'relaxation_velocity_g': r'气相速度松弛因子至([\d.]+)',
            'relaxation_velocity_l': r'液相速度松弛因子至([\d.]+)',
        }
        
        # 按迭代分割
        iter_blocks = re.split(r'第\s+\d+\s+次迭代:', content)[1:]  # 跳过第一个空块
        
        for idx, block in enumerate(iter_blocks, start=1):
            self.data['iterations'].append(idx)
            
            # 提取残差
            for key, pattern in patterns.items():
                if key == 'iteration':
                    continue
                
                # 处理对角占优度的特殊情况（可能有两个值：vapor和liquid）
                if key == 'diag_dominance':
                    matches = re.findall(pattern, block)
                    if matches:
                        # 取第一个值（通常是vapor的）
                        self.data['diag_dominance_vapor'].append(float(matches[0]))
                        if len(matches) > 1:
                            self.data['diag_dominance_liquid'].append(float(matches[1]))
                        else:
                            if len(self.data['diag_dominance_liquid']) > 0:
                                self.data['diag_dominance_liquid'].append(self.data['diag_dominance_liquid'][-1])
                            else:
                                self.data['diag_dominance_liquid'].append(np.nan)
                    else:
                        for dk in ['diag_dominance_vapor', 'diag_dominance_liquid']:
                            if len(self.data[dk]) > 0:
                                self.data[dk].append(self.data[dk][-1])
                            else:
                                self.data[dk].append(np.nan)
                    continue
                    
                match = re.search(pattern, block)
                if match:
                    value = float(match.group(1))
                    self.data[key].append(value)
                else:
                    # 如果没有匹配，使用上一个值或NaN
                    if len(self.data[key]) > 0:
                        self.data[key].append(self.data[key][-1])
                    else:
                        self.data[key].append(np.nan)
        
        print(f"解析完成！共{len(self.data['iterations'])}次迭代")
        self._create_dataframe()
        
    def _create_dataframe(self):
        """创建pandas DataFrame方便分析"""
        self.df = pd.DataFrame(self.data)
        
    def plot_residuals(self, save_path: str = None):
        """
        绘制残差变化图
        
        Args:
            save_path: 保存路径（可选）
        """
        fig, axes = plt.subplots(2, 2, figsize=(15, 10))
        fig.suptitle('收敛残差变化趋势', fontsize=16, fontweight='bold')
        
        # 1. 轴向动量残差
        ax1 = axes[0, 0]
        ax1.semilogy(self.data['iterations'], self.data['res_vapor_axial'], 
                     'b-o', label='气相轴向', markersize=3, linewidth=1.5)
        ax1.semilogy(self.data['iterations'], self.data['res_liquid_axial'], 
                     'r-s', label='液相轴向', markersize=3, linewidth=1.5)
        ax1.axhline(y=5e-3, color='b', linestyle='--', alpha=0.5, label='气相目标(5e-3)')
        ax1.axhline(y=1e-3, color='r', linestyle='--', alpha=0.5, label='液相目标(1e-3)')
        ax1.set_xlabel('迭代次数')
        ax1.set_ylabel('残差 (对数刻度)')
        ax1.set_title('轴向动量残差')
        ax1.legend()
        ax1.grid(True, alpha=0.3)
        
        # 2. 横向动量残差
        ax2 = axes[0, 1]
        ax2.semilogy(self.data['iterations'], self.data['res_vapor_lateral'], 
                     'b-o', label='气相横向', markersize=3, linewidth=1.5)
        ax2.semilogy(self.data['iterations'], self.data['res_liquid_lateral'], 
                     'r-s', label='液相横向', markersize=3, linewidth=1.5)
        ax2.axhline(y=1e-2, color='b', linestyle='--', alpha=0.5, label='气相目标(1e-2)')
        ax2.axhline(y=5e-3, color='r', linestyle='--', alpha=0.5, label='液相目标(5e-3)')
        ax2.set_xlabel('迭代次数')
        ax2.set_ylabel('残差 (对数刻度)')
        ax2.set_title('横向动量残差')
        ax2.legend()
        ax2.grid(True, alpha=0.3)
        
        # 3. 质量守恒残差
        ax3 = axes[1, 0]
        ax3.semilogy(self.data['iterations'], self.data['res_mass_vapor'], 
                     'b-o', label='气相质量', markersize=3, linewidth=1.5)
        ax3.semilogy(self.data['iterations'], self.data['res_mass_liquid'], 
                     'r-s', label='液相质量', markersize=3, linewidth=1.5)
        ax3.axhline(y=1e-3, color='g', linestyle='--', alpha=0.5, label='目标(1e-3)')
        ax3.set_xlabel('迭代次数')
        ax3.set_ylabel('残差 (对数刻度)')
        ax3.set_title('质量守恒残差')
        ax3.legend()
        ax3.grid(True, alpha=0.3)
        
        # 4. 压力和空泡份额残差
        ax4 = axes[1, 1]
        ax4.semilogy(self.data['iterations'], self.data['res_pressure'], 
                     'g-o', label='压力', markersize=3, linewidth=1.5)
        ax4.semilogy(self.data['iterations'], self.data['res_void_fraction'], 
                     'm-s', label='空泡份额', markersize=3, linewidth=1.5)
        ax4.axhline(y=1e-2, color='g', linestyle='--', alpha=0.5, label='压力目标(1e-2)')
        ax4.axhline(y=1e-4, color='m', linestyle='--', alpha=0.5, label='空泡目标(1e-4)')
        ax4.set_xlabel('迭代次数')
        ax4.set_ylabel('残差 (对数刻度)')
        ax4.set_title('压力和空泡份额残差')
        ax4.legend()
        ax4.grid(True, alpha=0.3)
        
        plt.tight_layout()
        
        if save_path:
            plt.savefig(save_path, dpi=300, bbox_inches='tight')
            print(f"残差图已保存到: {save_path}")
        
        plt.show()
        
    def plot_diagnostics(self, save_path: str = None):
        """
        绘制诊断信息图
        
        Args:
            save_path: 保存路径（可选）
        """
        fig, axes = plt.subplots(2, 2, figsize=(15, 10))
        fig.suptitle('求解器诊断信息', fontsize=16, fontweight='bold')
        
        # 1. 对角占优度
        ax1 = axes[0, 0]
        if len([x for x in self.data['diag_dominance_vapor'] if not np.isnan(x)]) > 0:
            ax1.plot(self.data['iterations'], self.data['diag_dominance_vapor'], 
                     'b-o', label='气相', markersize=3, linewidth=1.5)
        if len([x for x in self.data['diag_dominance_liquid'] if not np.isnan(x)]) > 0:
            ax1.plot(self.data['iterations'], self.data['diag_dominance_liquid'], 
                     'r-s', label='液相', markersize=3, linewidth=1.5)
        ax1.axhline(y=1, color='k', linestyle='--', alpha=0.5, label='对角占优临界值')
        ax1.set_xlabel('迭代次数')
        ax1.set_ylabel('对角占优度')
        ax1.set_title('矩阵对角占优度 (>1为对角占优)')
        ax1.legend()
        ax1.grid(True, alpha=0.3)
        ax1.set_yscale('log')
        
        # 2. 松弛因子变化 - 压力
        ax2 = axes[0, 1]
        if len([x for x in self.data['relaxation_pressure'] if not np.isnan(x)]) > 0:
            valid_idx = [i for i, x in enumerate(self.data['relaxation_pressure']) if not np.isnan(x)]
            valid_iter = [self.data['iterations'][i] for i in valid_idx]
            valid_vals = [self.data['relaxation_pressure'][i] for i in valid_idx]
            ax2.plot(valid_iter, valid_vals, 'g-o', markersize=5, linewidth=2)
        ax2.set_xlabel('迭代次数')
        ax2.set_ylabel('松弛因子')
        ax2.set_title('压力松弛因子调整')
        ax2.grid(True, alpha=0.3)
        ax2.set_ylim([0, 0.2])
        
        # 3. 松弛因子变化 - 气相速度
        ax3 = axes[1, 0]
        if len([x for x in self.data['relaxation_velocity_g'] if not np.isnan(x)]) > 0:
            valid_idx = [i for i, x in enumerate(self.data['relaxation_velocity_g']) if not np.isnan(x)]
            valid_iter = [self.data['iterations'][i] for i in valid_idx]
            valid_vals = [self.data['relaxation_velocity_g'][i] for i in valid_idx]
            ax3.plot(valid_iter, valid_vals, 'b-o', markersize=5, linewidth=2, label='气相')
        ax3.set_xlabel('迭代次数')
        ax3.set_ylabel('松弛因子')
        ax3.set_title('气相速度松弛因子调整')
        ax3.grid(True, alpha=0.3)
        ax3.set_ylim([0, 0.5])
        
        # 4. 松弛因子变化 - 液相速度
        ax4 = axes[1, 1]
        if len([x for x in self.data['relaxation_velocity_l'] if not np.isnan(x)]) > 0:
            valid_idx = [i for i, x in enumerate(self.data['relaxation_velocity_l']) if not np.isnan(x)]
            valid_iter = [self.data['iterations'][i] for i in valid_idx]
            valid_vals = [self.data['relaxation_velocity_l'][i] for i in valid_idx]
            ax4.plot(valid_iter, valid_vals, 'r-o', markersize=5, linewidth=2, label='液相')
        ax4.set_xlabel('迭代次数')
        ax4.set_ylabel('松弛因子')
        ax4.set_title('液相速度松弛因子调整')
        ax4.grid(True, alpha=0.3)
        ax4.set_ylim([0, 0.5])
        
        plt.tight_layout()
        
        if save_path:
            plt.savefig(save_path, dpi=300, bbox_inches='tight')
            print(f"诊断图已保存到: {save_path}")
        
        plt.show()
        
    def plot_convergence_summary(self, save_path: str = None):
        """
        绘制收敛总结图：所有残差在一张图上
        
        Args:
            save_path: 保存路径（可选）
        """
        plt.figure(figsize=(12, 8))
        
        plt.semilogy(self.data['iterations'], self.data['res_vapor_axial'], 
                     'b-', label='气相轴向动量', linewidth=2, alpha=0.8)
        plt.semilogy(self.data['iterations'], self.data['res_liquid_axial'], 
                     'r-', label='液相轴向动量', linewidth=2, alpha=0.8)
        plt.semilogy(self.data['iterations'], self.data['res_vapor_lateral'], 
                     'b--', label='气相横向动量', linewidth=1.5, alpha=0.6)
        plt.semilogy(self.data['iterations'], self.data['res_liquid_lateral'], 
                     'r--', label='液相横向动量', linewidth=1.5, alpha=0.6)
        plt.semilogy(self.data['iterations'], self.data['res_mass_vapor'], 
                     'c-', label='气相质量', linewidth=2, alpha=0.7)
        plt.semilogy(self.data['iterations'], self.data['res_mass_liquid'], 
                     'm-', label='液相质量', linewidth=2, alpha=0.7)
        plt.semilogy(self.data['iterations'], self.data['res_pressure'], 
                     'g-', label='压力', linewidth=2, alpha=0.7)
        plt.semilogy(self.data['iterations'], self.data['res_void_fraction'], 
                     'orange', label='空泡份额', linewidth=2, alpha=0.7)
        
        # 添加收敛目标线
        plt.axhline(y=1e-3, color='gray', linestyle=':', alpha=0.5, label='1e-3')
        plt.axhline(y=1e-2, color='gray', linestyle='-.', alpha=0.5, label='1e-2')
        
        plt.xlabel('迭代次数', fontsize=12)
        plt.ylabel('残差 (对数刻度)', fontsize=12)
        plt.title('收敛过程总览 - 所有残差', fontsize=14, fontweight='bold')
        plt.legend(loc='best', fontsize=10, ncol=2)
        plt.grid(True, alpha=0.3, which='both')
        
        plt.tight_layout()
        
        if save_path:
            plt.savefig(save_path, dpi=300, bbox_inches='tight')
            print(f"收敛总结图已保存到: {save_path}")
        
        plt.show()
        
    def print_summary(self):
        """打印收敛统计摘要"""
        print("\n" + "="*60)
        print(" " * 20 + "收敛分析摘要")
        print("="*60)
        
        print(f"\n总迭代次数: {len(self.data['iterations'])}")
        
        # 最终残差
        print("\n最终残差值:")
        print(f"  气相轴向动量: {self.data['res_vapor_axial'][-1]:.3e} (目标: < 5e-3)")
        print(f"  液相轴向动量: {self.data['res_liquid_axial'][-1]:.3e} (目标: < 1e-3)")
        print(f"  气相横向动量: {self.data['res_vapor_lateral'][-1]:.3e} (目标: < 1e-2)")
        print(f"  液相横向动量: {self.data['res_liquid_lateral'][-1]:.3e} (目标: < 5e-3)")
        print(f"  气相质量守恒: {self.data['res_mass_vapor'][-1]:.3e} (目标: < 1e-3)")
        print(f"  液相质量守恒: {self.data['res_mass_liquid'][-1]:.3e} (目标: < 1e-3)")
        print(f"  压力:         {self.data['res_pressure'][-1]:.3e} (目标: < 1e-2)")
        print(f"  空泡份额:     {self.data['res_void_fraction'][-1]:.3e} (目标: < 1e-4)")
        
        # 收敛检查
        print("\n收敛检查:")
        converged = {
            '气相轴向动量': self.data['res_vapor_axial'][-1] < 5e-3,
            '液相轴向动量': self.data['res_liquid_axial'][-1] < 1e-3,
            '气相横向动量': self.data['res_vapor_lateral'][-1] < 1e-2,
            '液相横向动量': self.data['res_liquid_lateral'][-1] < 5e-3,
            '气相质量': self.data['res_mass_vapor'][-1] < 1e-3,
            '液相质量': self.data['res_mass_liquid'][-1] < 1e-3,
            '压力': self.data['res_pressure'][-1] < 1e-2,
            '空泡份额': self.data['res_void_fraction'][-1] < 1e-4,
        }
        
        for name, status in converged.items():
            symbol = "✓" if status else "✗"
            print(f"  [{symbol}] {name}")
        
        all_converged = all(converged.values())
        print(f"\n整体收敛状态: {'✓ 已收敛' if all_converged else '✗ 未完全收敛'}")
        print("="*60 + "\n")


def main():
    """主函数：演示如何使用分析器"""
    print("=" * 60)
    print(" " * 15 + "收敛分析快速工具")
    print("=" * 60 + "\n")
    
    # 查找最新日志
    log_file = find_latest_log()
    
    if log_file is None:
        print("❌ 未找到日志文件！")
        print("请确保在以下目录中有日志文件：")
        print("  - e:\\mycode\\subchannel\\TwoFluidSCA\\results")
        print("  - e:\\mycode\\subchannel\\TwoFluidSCA\\src\\results")
        sys.exit(1)
    
    print(f"📁 找到最新日志: {log_file.name}")
    print(f"📅 修改时间: {Path(log_file).stat().st_mtime}")
    print()
    
    # 创建分析器
    analyzer = ConvergenceAnalyzer(str(log_file))
    
    # 解析日志
    analyzer.parse_log()
    
    # 打印摘要
    analyzer.print_summary()
    
    # 生成图表
    results_dir = log_file.parent / "analysis_plots"
    results_dir.mkdir(exist_ok=True)
    
    print("\n📊 正在生成可视化图表...")
    
    try:
        analyzer.plot_convergence_summary(str(results_dir / "convergence_summary.png"))
        print(f"  ✓ 收敛总结图: {results_dir / 'convergence_summary.png'}")
    except Exception as e:
        print(f"  ✗ 收敛总结图生成失败: {e}")
    
    try:
        analyzer.plot_residuals(str(results_dir / "residuals.png"))
        print(f"  ✓ 残差分析图: {results_dir / 'residuals.png'}")
    except Exception as e:
        print(f"  ✗ 残差分析图生成失败: {e}")
    
    try:
        analyzer.plot_diagnostics(str(results_dir / "diagnostics.png"))
        print(f"  ✓ 诊断信息图: {results_dir / 'diagnostics.png'}")
    except Exception as e:
        print(f"  ✗ 诊断信息图生成失败: {e}")
    
    print(f"\n✅ 所有图表已保存到: {results_dir}")
    print("\n分析完成！")


if __name__ == "__main__":
    main()
